Quantum Statistical Physics of Glasses at Low Temperatures 



o 
o 

o 

Q 

On 



J. van Baardewijk and R. Kvihn. 

Department of Mathematics, King's College University-London, Strand, London, WC2R.2LS, UK. 

(Dated: December 19, 2009) 

Wc present a quantum statistical analysis of a microscopic mean-field model of structural glasses at low temper- 
atures. The model can be thought of as arising from a random Born von Karman expansion of the full interaction 
potential. The problem is reduced to a single-site theory formulated in terms of an imaginary-time path integral 
using replicas to deal with the disorder. We study the physical properties of the system in thermodynamic equi- 
librium and develop both pcrturbative and non-perturbativc methods to solve the model. The perturbation theory 
is formulated as a loop expansion in terms of two-particle irreducible diagrams, and is carried to three-loop order 
in the effective action. The non-perturbativc description is investigated in two ways, (i) using a static approxi- 
mation, and (ii) via Quantum Monte Carlo simulations. Results for the Matsubara correlations at two-loop order 
perturbation theory are in good agreement with those of the Quantum Monte Carlo simulations. Characteristic 
low-temperature anomalies of the specific heat arc reproduced, both in the non-perturbativc static approximation, 
and from a three-loop pcrturbative evaluation of the free energy. In the latter case the result so far relies on 
using Matsubara correlations at two-loop order in the three-loop expressions for the free energy, as self-consistent 
Matsubara correlations at three-loop order are still unavailable. Wc propose to justify this by the good agreement 
of two-loop Matsubara correlations with those obtained non-pcrturbatively via Quantum Monte Carlo simulations. 
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I. INTRODUCTION 

Glasses are known to exhibit distinctive low- 
temperature properties that differ substantially from 
those of crystalline solids and are referred to as glassy 
low-temperature anomalies. For instance, at low tem- 
peratures the specific heat and thermal conductivity in 
crystals show a familiar Independence. In glasses the 
specific heat is found to increase approximately linearly 
with the temperature at T < IK, while the thermal 
conductivity increases approximately as T 2 in this low 
temperature range [l[. At higher temperatures between 
1 and 20 K the thermal conductivity is approximately 
constant, while the specific heat C shows a peak when 
displayed as C/T 3 , usually referred to as the Bose-peak. 
In phenomenological models such as the Standard Tun- 
nelling Model [IjH and the soft potential model [H1H one 
postulates that a broad spectrum of tunnelling centres is 
responsible for the properties at T < 1 K . 

Surprisingly, the glassy anomalies show a noticeable 
degree of universality at T < IK, whereas between ap- 
proximately 1 and 20 K these depend more on the spe- 
cific materials [l|. There are currently two main con- 
tending theories to explain this fact. Following ideas 
of Yu and Leggett @, it has been suggested as result- 
ing from a collective effect due to interactions between 
the tunnelling excitations 0] ■ Alternatively, it is thought 
to be a property of the potential energy landscape cre- 
ated by glassy freezing at high temperatures. This also 
defines a phenomenon of collective origin but involves 
no quantum effects Q. Universality in the second in- 
terpretation is understood as a result of separation of 
the energy scales involved in glassy freezing on the one 
hand side, and those relevant for the low-temperature 
phenomena on the other hand side [8, 9]. Whereas the 
existence of tunnelling centres is part of the initial as- 
sumptions in [3], these are shown to arise naturally as 
a result of microscopic interactions in the model glasses 



studied in @, Q , and do indeed give rise to the charac- 
teristic low-temperature anomalies. The work in Q is 
perhaps appropriately characterised as a strong coupling 
approach to glassy low-temperature physics. A comple- 
mentary weak coupling approach [13|, [lj] to the same 
phenomena takes weak residual interactions between a 
set of quasi-local collective modes as starting point and 
describes low-temperature anomalies in terms of a vibra- 
tional instability occuring in systems of this type. 

The analysis in [H, [t| is still semi-classical in the sense 
that it is based on an analysis of quantum effects in a 
glassy potential energy landscape whose properties were 
determined via classical statistical mechanics. The aim 
of the present paper is to overcome this deficit and study 
the system in a full quantum statistical formulation right 
from the outset. Focus will be here on the translationally 
invariant model proposed in [9j. 

We shall proceed along the lines of general methods de- 
veloped for quantum spin-glasses. In particular, we apply 
the Matsubara formalism to construct an imaginary-time 
path integral representation of the partition function and 
the replica-method to deal with the disorder. The sites 
are decoupled by introducing order parameters for which 
the functional integral is evaluated by the method of 
steepest descent. The result is an effective single-site 
theory and a set of functional self-consistency relations 
for the order parameters. These methods are similar to 
those used for models studying spin-glass transitions in 
quantum spin-glasses. Examples are the SK-model of 
spin-glasses generalised to quantum s pins II Oil and the 
quantum spherical p-spin glass model [lj, OjJ . Here we 
shall not concern ourselves with the glass transition but 
concentrate on evaluating the physical properties at low 
temperatures and in particular the specific heat anomaly 
in the 1 K region. 

To solve the effective single-site theory we first apply a 
perturbative method in terms of two-particle irreducible 
(2PI) diagrams, which is based on an expansion in pow- 
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ers of the full interacting correlation functions. This 
amounts to summing infinite classes of diagrams and can 
therefore also capture effects of a non-perturbative na- 
ture. After this we develop a non-perturbative theory 
proper. The result is a set of functional self-consistency 
equations for the order parameters which we first treat 
with Quantum Monte Carlo simulations. 

Following this we construct a solvable version of the 
non-perturbative theory, using a simple approximation 
known from quantum spin-glass theory as the static ap- 
proximation. This scheme was first introduced as a vari- 
ational Ansatz in (T(jj where the time-dependent order 
parameter was approximated by a time-independent con- 
stant. 

Given the complications of dealing with quantum 
fluctuations in this model, we presently restrict the 
analysis of both perturbative and non-perturbative 
theory to the replica symmetric approximation. In 
support of this we mention that the effects of replica 
symmetry breaking on the low-temperature anomalies 
were found to be small at the semi-classical level [1, []| . 

This paper is organised as follows. In section |TT] we 
summarise the main ingredients of the proposed glass 
model. In section ITO1 we give the many-particle partition 
function represented by an imaginary-time path integral 
and introduce replicas to handle the disorder averaging. 
Section HVl gives an account of the effective single-site for- 
mulation, deriving the effective action and the functional 
self-consistency relations for the order parameters. Then 
in section IVl we treat perturbative and in section [VI] non- 
perturbative solution methods. The numerical results for 
the order parameters and the specific heat are discussed 
in section Ivm Finally, the conclusions are drawn in sec- 
tion EH 



have elements of randomness and frustration. There are 
two possibilities to make analytic progress: use a mean- 
field approximation for a microscopically 'semi-realistic' 
model or alternatively formulate a model for which such 
a mean-field approximation would be exact. The latter is 
the approach we have taken here, a justification of which 
we believe is provided by the results. The absence of 
phonons is of course one of the unavoidable consequences 
of adopting a mean-field approximation. 

Following *9|, the potential function is taken to repre- 
sent the first terms of a Born von Karman expansion of 
the full interaction energy about the reference positions. 
A further requirement of global ^-symmetry (u <-> — u) 
excludes the odd orders in this expansion. Matters are 
further simplified by taking the Ui to be scalar, resulting 
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Jij (Ui - Uj) 2 + — (Ui - Uj) 



(2) 



The glassy properties are represented at the quadratic 
level in @, defining a random-interaction term with ran- 
dom interaction strengths Jij. The quartic term (taken 
to be non-random) is necessary in order to stabilise the 
system as a whole, and so g > 0. The parameters 
are quenched and taken independent with equal Gaus- 
sian distribution Af(0, J 2 /N) for each combination 
The 1/N scaling of the variance and the quartic interac- 
tion term in ^ ensures that the thermodynamic energy 
is proportional to N. The construction as presented here 
allows the system to be analysed within replica mean- 
field theory, similar to that of the SK-model for spin- 
glasses [Hi, its generelization to quantum spin-glasses 
[lOj ] and quantum spherical p-spin glasses (Til Il2| . 



II. THE GLASS MODEL 

Starting point is the microscopic model for a glass- 
like system at low temperatures, proposed in Q. To 
summarise its main ingredients, consider a system of N 
degrees of freedom (called particles) with the following 
model-Hamiltonian 



JV 2 
2=1 



(1) 



where the pi denote the momenta and the m, the masses, 
which for simplicity are taken equal for each particle. The 
variables u = (u\, ..,ttjv) represent coordinate deviations 
from pre- assigned reference positions. Since glassy low- 
temperature physics is universal, the interaction poten- 
tial V(u) need not contain details of the specific atoms 
and their specific interactions and is taken as simple as 
possible, yet containing enough detail to reproduce glassy 
low-temperature physics. The minimum requirement is 
that it should respect global translation invariance and 



III. THE PARTITION FUNCTION 

The quantum statistical partition function for the fixed 
disorder configuration {Jij} in a basis of coordinate 
states |ti) = \ui) . . . \uf?) is 

Zj = Tr exp(-/3H) = [ du (u\cxp(-f3H)\u), (3) 



where the Hamiltonian H is defined by ([I} with itj and 
Pi replaced by the operators Ui and pi . In the Matsubara 
formalism the path integral representation of Q is con- 
structed using the Lie- Trotter product formula [3, [TtJ 



exp(-/3T - f3V) 



lim 

r — >oo 



{ ex p(- v) exp ( 



(4) 

Insertions of 1 = / du k \u k ){u k \ and 1 = J dp k \p k )(p k \, 
together with definitions of imaginary time r k = kAr 
(k = 0,..,r — 1) and time-step At = leads to the 
following path integral representation of (J3]) 



Zj = / Du exp 



" H A[U] 



(•5) 
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where the integration is in the functional sense with a 
measure defined as 



r-l N 



fc=0 i=l 



2irh 2 P 



dUi(T k ). 



(0) 



The functions u,(r) satisfy the periodicity conditions 
Ui(Q) = Ui(h(3). The Euclidean action reads 



A[i 



Tui 



m (dui(r)\ 2 



dr 



+ V(u(t)) . (7) 



The interaction potential V(u(t)) equals the expression 
in |(3J) with Ui replaced by Ui(r). 

In order to study the equilibrium properties of the 
model we need to compute the disorder averaged free 
energy density /. The replica trick [l8| allows us to eval- 
uate this as 



1 



1 



(if = -log Z., = lim — log (Z, ; y\ (8) 
W n->o Nn 



where the overline denotes the average over all realiza- 
tions of the random interaction. To end this section we 
list the expression for (Zj) n , i.e. the replicated version 
of© 



/n n 
l[Vu a exp(--]T.4[« a ] 
,-, — 1 — 1 



(9) 



The index a numbers the replicas and u a 



1 ' ") u Nl 



IV. EFFECTIVE SINGLE-SITE FORMULATION 



In order to evaluate ((8]) we first average over all re- 
alizations of the random potential. This is achieved by 
carrying out the independent Gaussian integrations over 
the set { Jij}- The result is 



( Zj )n = /n^cxp[^{^(^/dr[<(r)-^(r)] 2 ) 2 - ^e/^km-^)] 4 } 

J a=l i<j a J a J 

x ex p[-^E/ dr f <( r ) : 



(10) 



The expansions of the powers in the first line of (|10[) con- This means that we do not have to consider terms of the 

tain many terms that vanish due to the following Ansatz. kind . ^2 ab J dr J dr' u^(t) 2 u^(t')u^(t'), which would 

We assume the global Z 2 -symmetry to remain unbroken complicate the formulation considerably. The result after 

after quantization: 1/iVJ^ a"(r) = for all r £ [0, ft/3]. the expansions is 



J exp[ W drdr'{ 

ia ab 



J 



J 2 



4/ , v vE<(^(-')) 2 + ^E<(^(-') 2 } 

i i 

2 + 3<w 4 ) + |(E<w 2 ) 2 }]- 

i 

i 



(ii) 



The sites are decoupled with two sets of Gaussian trans- where Z e s defines the effective single-site partition func- 
formations after which becomes tion 



(ZjY = J V{q aa (T,T)}V{q ab (T,T')} 

exp{N(-^X[q] + logZ off )}, (12) 



z ° a = J n Vua exp ( ~ \ Scs w 



(13) 
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The non-fluctuating part in (|12|) is denned as 

-3gJ2 f dTq aa (r,T) 2 . (14) 

The effective single-site action reads 

5 c ff[u Q ] = ^^ y rfrdr' u a {T)q^ ab (T, t')«6(t') +5 in t[u ]- 

The interaction part iS; nt [u a ] contains a quartic term non- 
local in time and quartic term local in time 

J 2 f 

Si n t[u a ] = -3-^^ / drdr' u a (T) 2 u b (T') 2 

ab 

+ .9^ fdru a (r)\ (16) 

The 'free' inverse propagator is 
d 2 

9o,ab( T > T ') = {~m— + 12 gq aa (T,T)}S ab S(t-t') 



dr 
J 2 

-T-gab(T, T 1 ) 



(17) 



Observe here that S e s [u a ] in (fl~3| and thus also Z e g, de- 
pend functionally on g a (,(r, t'). 

The replicated partition function (fT2|) can be treated 
with the saddle-point method. At the saddle-points we 
have 

(^~exp{AT(- ±X[q] + logZ cff )}, (18) 

where the saddle-point fields q a b( T ,T') are the order pa- 
rameters of the theory. The saddle-point equations result 
in the following functional self-consistency relations for 
the order parameters 



q a b(T,T') = (u„(r) u b (r')) 



(19) 



The angular brackets (...) denote the quantum thermo- 
dynamical average mediated by the effective action (fl~5|) . 
The order parameters q a b(T, r') are the full interacting 
correlation functions of the single-site theory, from here 
on called Matsubara correlations. 

The Matsubara correlations are time-translational in- 
variant since we are studying an equilibrium problem. 
They are also symmetric in time due to the time-reversal 
invariance of the action (1151). i.e. we have 



q ab {r, t') = q ab (T - t') = q ab {r' - r). (20) 

Furthermore the Matsubara correlations q a b{r — t ) are 
HP time-periodic. 

We should mention that the first interaction term in 
(I16[) defines a complete square, which could be linearised 



at the cost of introducing a Gaussian family of systems. 
However, we have chosen not to do this at this stage. It 
would lead to more complicated saddle-point equations 
for the order parameters when solving the single-site the- 
ory perturbativcly. However, we shall linearise this inter- 
action term in the non-perturbative treatment. 



V. PERTURBATION THEORY 



2PI-effective action formalism 



To solve the single-site theory perturbatively, we need 
a formalism that expands the path integral Q13[) in terms 
of a further effective (classical) action r c ff[q] 



Z e s[q] = j Y[ Vu a, exp( - ^S cS [u a .q] 

>(-~r eff [g]). (21) 



exp 



Here we have explicitly referred to the functional depen- 
dences on the Matsubara correlations in S e f[[u a ,q] and 
in Z e ff[g]. Remember, this dependence is due to the ap- 
pearance of q a b(j, t') in the inverse propagator (|17[) . The 
qab(T, t') define the full interacting correlations as we saw 
in the previous section. As regards to a perturbative ex- 
pansion of the path integral, the most efficient way is to 
also express r e ff[q] entirely in terms of the full interacting 
correlations, which was already assumed in the notation 
in (f2"Tj) . For this we choose the two-particle irreducible 
(2PI) effective action approach, developed in field theory 
[19L |20|. which is indeed based on an expansion of r e g[g] 
in powers of the full interacting correlators q a b(T, t') and 
involves only 2PI diagrams. The 2PI nature of the di- 
agrams has the additional advantage of considerably re- 
ducing the number of diagrams that need to be included 
in the expansion. Also, as the expansion is in terms of 
the full interacting correlators, the 2PI approach effec- 
tively amounts to summing infinite classes of diagrams 
of a conventional perturbation expansion, thus enabling 
it to capture effects of a non-perturbative nature. Use- 
fulness of the 2PI effective action approach for the s tudy 
of stochastic dynamical systems was advocated in [2l| . 
Its application to the analysis of glassy systems was sug- 
gested in (22 |. 

Before presenting the series expansion of r o ff[<7], we 
discuss the key ingredients of the 2PI effective action 
approach, in a formulation appropriate for the present 
problem. Following [19(, one first adds a two-body source 
term to the action S c s- This defines a generating func- 
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tional 

Zes[K, q\ = I Y\. ' Dua exp ( - \ { ^ cff t Ua ' ^ 



dr dr 



1 u a (T)K ab (T, r')ub(r') j ) 
exp(-iw eff [X,?]) ! (22) 



giving 



SW cS [K, q] 



1 



<5# q6 (t,t') 2 
From and dH]) follows 



Sqab(T,T') hSq ab (T,T') 



= 0, 



(23) 



(24) 



We shall need this in the equations of motion for 
qab( T i T ')i to be constructed next. In order to eliminate 
K in favour of the full interacting correlations q, one per- 
forms the following Legendre transformation 

TeffM = W cS [K,q] drdT'q ab (T,T')K ab (T,r'), 

ab 

(25) 



giving 



dT eS [q} _5W cS [K,q] 1 



Sq a b(r,T') Sq ab (T,T') 



K ab (r, t'). (26) 



Then setting the source field K ab (r, r') to zero and us- 
ing ([24]) results in the following 'equations of motion' for 

QabiT, T') 



ST ctf [q} 1 6X[q] 



5q a b(T,T') hSq ab (T,T') 



0. 



(27) 



Next is to describe the series expansion of r R ffJql, which 
we present in the standard form as derived in [l9l . l20l | 

T cff [G] = ^Trqo 1 G+ ^TrlogG- 1 + £r p [G], (28) 

with Green's function G ab (T,r') = q ab (T,T')/h. As al- 
ready mentioned, the variable (fab ( T > r ') defines the full 
interacting correlation functions and satisfies the equa- 
tions of motion (|27|) . The traces are taken in the func- 
tional sense. The first two terms in (|28|) define what is 
called the one-loop contribution. The terms denoted by 
T p [G] define the p-loop contributions. These are repre- 
sented by 2PI diagrams containing p loops. In the next 
section we shall discuss the rules for constructing such 
diagrams. A diagram is said to be 2PI if it does not 
become disconnected upon cutting two lines. The fact 
that the contributions J2^L2^p[G] define 2PI diagrams 



is understood from the following argument. The equa- 
tions of motion (|27| with substitutions of (fl"4|) and (|28|) 
just result in the Dyson equations 



G~ 



^(EW + ^-iG] 



i 



SG\ 



(29) 



Since the second part of (|2"9"|) defines the proper self- 
energy, the diagrams of which are known to be one- 
particle-irreducible, clearly the diagrams of the terms 
r p [G] must be 2PI, (the diagrams of X[G] are 2PI). In 
the expansion of 2PI diagrams we shall consider two or- 
ders, the two-loop and the three-loop order. They are 
analysed in section IV CI and TV Dl below. 

One should mention that (|2"8|) is only valid for the spe- 
cial case (u a (r)) = 0. If (w q (t)) = U a ^ a further 
source term in (|22|) is needed. This would lead to ex- 
tra terms depending on U a in the effective action (|28l) . 
and a further equation of motion ST e ff/SU a = [19j . 
We shall not concern ourselves with this case since we 
have assumed global symmetry to remain unbroken, i.e. 
(u a (r)) = (see section ITV)) . 

In the 2PI effective action formalism the classical form 
(|2~Tj) results in the following expression for free energy 



1 



lim l — (X[q]+T oS [q] 
n—>o nn 



where we have used dH]) and (|T8l) . 



(30) 



B. Rules for 2PI diagrams 



The two interaction terms in (|16p determine the fol- 
lowing rules for constructing 2PI diagrams and their ex- 
pressions. Vertices are labelled by a replica index and 
a time variable. Two vertices labelled by (a, r) and 
(b, t 1 ) are connected by a solid line contributing a fac- 
tor HG ab (T, t') — q a b{T, t'). The two interaction terms 
define two types of vertices. Firstly, there is a g-type 
vertex which is represented by a dot • and contributes 
a factor —g/h. Secondly, there is a J-type vertex which 
is represented by a cross x . The J-type vertices always 
come in sets of two. The two are connected by a dashed 

line x x . A dashed line contributes a factor J 2 /8h 2 . 

Each diagram gets an extra factor —h due to the prefac- 
tor —1/h in (|13j) . Then there is a further permutational 
factor to consider which we write in front of the diagram. 
Next is to collect all factors of a diagram and multiply 
them. Finally, one needs to sum over all replica indices 
and integrate over all time variables. 

When constructing diagrams of order p, all distinct di- 
agrams containing p loops are added together, resulting 
in the final expression for r p [g]. When counting loops, 
both solid lines and dashed lines need to be considered. 
The number of loops in a diagram is equal to the power 
of H in its expression. A subtlety here is that when count- 
ing powers of h, one factor h in the contribution of each 
dashed line J 2 /8h 2 must not be taken into account. This 
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is the factor h that originates from the coupling constant 
J 2 /8h of the non-local interaction term in (fTo| . which 
represents an interaction parameter as a whole. Finally, 
it should be noted that each expansion order T p [q] con- 
tains diagrams of 0(n) in replica and diagrams of 0(n 2 ) 
or higher order. Since n — > the latter can be ignored. 
The 2PI diagrams for two-loop order (p = 2) and three- 
loop order (p — 3) are given in Fig. [TJ The prefactors 
here result from permutations. 




FIG. 1: Two-loop and three-loop 2PI diagrams 



Again, we derive the equations for the three-loop Matsub- 
ara correlations q a b{ T T T ') from the equations of motion 
(|27| . resulting in 

l h J 2 

+6g5 ab 6(T - t') q aa (T, t) - ^j—q a b{T, t') 3 
+ ^V(t, r 1 ) J dr" 9ac (r, r"f . (34) 

Solutions of the self-consistency equations (|3"2"|) and 
(I34p were attempted in a replica symmetric approxima- 
tion, and using Fourier transform techniques as described 
below. First, however let us turn to a description of the 
two non-perturbative approaches we have looked at, both 
also in a replica symmetric version. 

VI. NON-PERTURBATIVE THEORY 



A. Replica symmetry 



C. Two-loop order 

The two-loop 2PI diagrams are shown in the first line 
of Fig. [1] According to the rules listed in section IV Bl 
these diagrams represent the following expressions 



J 2 f 

?2[q] = -Jkz2 drdT ' qab(r, 

ab 

+ 3.9^ J dr q aa (r, rf 



J\2 



(31) 



where we have used g b(r, t') = %G ab (T,T'). The equa- 
tions for the two-loop Matsubara correlations q ab {T,T') 
are derived from the equations of motion (|27p . resulting 
in 

l h J 2 

= 2^afc( T > T ')- 2 q ab( T ' T ')- ^1^(t,t') 

+6gS ab S(T-T')q aa (T,T). (32) 

These were solved in a replica symmetric approximation, 
and using Fourier transform techniques as described in 
section IVIIB II below. 



D. Three-loop order 

The three-loop 2PI diagrams are shown in the second 
line of Fig. [T] These diagrams represent the following 
expressions 

12 g 2 



r 3 M 



3J 2 g 

h 2 



£ J drdr' q ab (T,r') 4 

ab 

J2 I dTdT'dr" qab (T,T') 2 q ac (r,T") 2 (33) 



The non-perturbative theory is constructed in the 
replica symmetric (RS) approximation. We mention here 
that the effects of replica symmetry breaking on the low- 
temperature anomalies have been small for the semi- 
classical treatments [8|, |9( . The RS form of the Matsubara 
correlations is 



qaa(r,T') = q d (T,T'), 

q a b(r,T') = q (for a ^ 6), 



(35a) 
(35b) 



independent of the replicas a, b. The off-diagonal Mat- 
subara correlation (|35b[) is assumed time-independent. 
The argument is that the replicas are independent and 
time-translationally invariant, and that the origin of time 
could be chosen independently for each replica [loj . 

The RS Ansatz allows us to decouple the replicas in the 
effective single-site action (fl~5|) . This concerns 2 terms. 
The first term to decouple is the non-local quartic inter- 
action given in (fTo| . We linearise this with a Gaussian 
variable z. The second term to consider comes from the 
part of qQ 1 in (fT7|) that involves q ab (T,T r ). This is de- 
coupled by means of the Gaussian variable z. After de- 
coupling, the effective single-site partition function (fT3]) 
becomes 

Z c s= [vzVz[[vuexp(^-^S R s[u-,z,z]^ , (36) 

where T>z — dzj\f^/K exp (— \z 2 ) and the same form for 
z. The decoupled RS action reads 



S RS [u;z,z] 



dr 



2 V dr 



- (l2g q d (T, r) + Jz^u(t) 2 + g u{rf 
(t,/) - ? )u(r)u(r')l (37) 



4 1 
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Substituting (|36| in (fT8| gives for the free energy |8]) 
PS = - f VzVz log{ /" 2? U exp( - i5 RS ) } 



with <7d(0) = gd(r, r). The free energy (|38|) becomes 



(38) 



where the non-fluctuating part is now defined as 

*RsM = ^1 drdr' (q d (r, r'f -q 2 ) -3 .g / o!t g d (r, r) 2 . 

(39) 

Finally, the RS saddle-point equations (11911 become 

q d { T ,T') ={{u{t)u{t'))) z - z , (40a) 



d =((u(t)) 2 ) z 



(40b) 



where (...) denotes an average mediated by the action 
iSrsj and (•••)«,£ the averages w.r.t. the Gaussians z and 
z. The theory defined by (|37M0p is non-local in Mat- 
subara time. There are currently no analytic techniques 
available to solve the self-consistency problem for the 
Matsubara correlations in this theory non-perturbatively 
while keeping the full complexity arising from this fact. 
Two different methods will be used to deal with this prob- 
lem, a numerical approach using Quantum Monte Carlo 
simulations as proposed in (23|, and the so-called static 
approximation due to Bray and Moore [icj |. 



B. Static approximation 

In this approximation the non-local term (the third 
line in (|3Tj) ) is approximated by taking for the diagonal 
Matsubara correlations a simple time-independent trial 
function 



Qd(T,r') = q d , 



(41) 



which may but need not assumed to be equal to qd(i~, r) = 
q d (0). Such static approximation scheme was introduced 
as a variational Ansatz in ^fj| when studying the spin- 
glass transition for the SK-model generalised to quantum 
spins. After applying (14111 we are able to linearise the 
non-local term of ([37]) by means of a further Gaussian 
transformation, defined by a Gaussian variable v. This 
results in the following static action 



5 st \u,v:z,z] 



dr 



m / du{r) \ 2 1 2 
2 V dr ) 2 V 
— JVC v u(t) + d\{z)u{r) 
+ d 2 {z)u{r) 2 + <?w(t) 4 1, (42) 



with C = P{qd — q). The random parameters d\ and d 2 
are defined as 



di(z) = -Jz^/q 

d 2 (z) = l(l2gq d (0) + Jz), 



VzVz\o%\ \ Vu 



dv 



exp 



1 



S*0J 

(44) 



(43) 



with the non-fluctuating part Ars [l] given by (|39p . 

Before presenting the self-consistency relations for the 
order parameters, let us first take a closer look at the 
static action (|42|) . The parameters di(z), d 2 (z) and g 
are the coupling constants of a potential energy for the 
system defined by the local quantum variable u(t). The 
(quenched) Gaussians z and z imply that we have here a 
heterogeneous family of such systems. This results in a 
broad spectrum of tunnelling and vibrational excitations, 
which we shall discuss in a separate section below. 

Interestingly, there is another term in the action 
namely —J\fC vu(t), which is of a different nature. The 
constant Jy/C defines a coupling constant for the bilin- 
ear interaction between the variable u(t) and the classi- 
cal (annealed) degree of freedom v. The latter is indeed 
a classical variable since it has no kinetic term associ- 
ated with it. The appearance of this 'annealed' degree 
of freedom is a result of formal mathematical analysis 
(as are the quenched variables z and z). Both are 'in- 
terpreted' and they acquire different meanings due to 
the different ways in which they appear in the theory. 
The z and z variables are frozen, resulting in an ensem- 
ble of double-well and single-well potentials, whereas the 
variable v defines a dynamical degree of freedom. The 
coupling to v is entirely analogous to the coupling to 
a heat-bath of phonons, as it is postulated in the phe- 
nomenological models of glassy low-temperature anoma- 
lies 0,0,13,0], though the details are of course different. 
Whereas phenomenological models postulate a coupling 
of local degrees of freedom to the strain-field of a heat- 
bath of phonons as an additional ingredient, the coupling 
to a harmonic classical variable v in the present case 
emerges through the (approximate) mathematical treat- 
ment of quantum fluctuations. Incidentally the presence 
of a heat-bath like background system could have been 
inferred directly from the appearance of retarded inter- 
actions in (|15p - (fT7| as such retarded interactions are the 
usual hallmark of effective descriptions of systems em- 
bedded into larger systems, after intergrating out the de- 
grees of freedom of those larger systems. 

Next we evaluate the free energy as a variational esti- 
mate w.r.t. the static Matsubara correlations qd(0),qd 
and q. The static formulation (|42H44[) defines a the- 
ory local in time. We were able to find numerical so- 
lutions for the Matsubara correlations in two different 
approaches, a three-variable approach in terms of the 
variables 9d(r, r) = qd(0),Qd and q, and a two-variable 
approach in terms of the variables qd and q, assuming 
9d(Q) — Qd- The variational equations in the three- 
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variable approach are 

q d (0) = ((u(r) 2 )) zz , 
rC ={{vu(r))) z - z , 
Q ={{u{t)) 2 ) z ,. 



J 



(45a) 
(45b) 
(45c) 



Here (...) denotes an average mediated by the static ac- 
tion JH]), while (...) z ,z denotes the Gaussian averages 
over z and z. In the two-variable approach the varia- 
tional equations are 

(ffc = 12 | ( {{u{rf)) z , -q d ) + J^( z { u (r))) zz , 

(46a) 

q = ((u(r)) 2 ) z - z . (46b) 



We will solve the functional self-consistency equations 
(|45|) and (|46]) reverting to an operator description, and 
using truncated Hilbert-spaces as described in section 
IVIIB 31 below. 



Tunnelling and vibrational excitations 



The potential energy in the action (|42| contains an 
ensemble of single- well and (asymmetric) double-well po- 
tentials because of the stochastic nature of the parame- 
ters di(z) and d,2(z) as defined in (|43|) . These single- 
well and double-well potentials are responsible for re- 
spectively vibrational excitations and the characteristic 
tunnelling excitations in the system [9j . We see that this 
potential structure arises naturally as a result of micro- 
scopic interactions defined by the model. In this context 
we mention the phenomcnological soft potential model 
[3, [B[ in which one postulates the existence of an ensemble 
of classical potentials V(u) = diu + d,2U 2 +gu 4 , providing 
a semi-classical analysis of its tunnelling- and vibrational 
states. Whereas that model assumes a uniform distribu- 
tion of the parameters d± and c?2 , the quantum statistical 
treatment of the microscopic model as presented here pre- 
dicts a Gaussian distribution of the parameters d\(z) and 
d,2(z). Furthermore, d\(z) and ££2(2) are parameterised 
by the disorder strength J and the order parameters q 
and <?(i(0). The collective nature of the latter can be 
seen as the origin of universality of the low-temperature 
physics predicted by the model. 



VII. NUMERICAL RESULTS 
A. Scaling 

For numerics and representation of results we represent 
the theory constructed in the previous sections in terms 
of dimensionless variables and parameters. Starting from 



a microscopic length scale uq we define the following en- 
ergy scales 



En = 



4 

gu 0l 



E, 



Ju , 



(47) 



where En defines the quantum energy scale. The dimen- 
sionless ratios of the variables are 



u _ q 
— ) 1 = — ■ 



Those for the parameters arc defined as 



E : 



J = 



Ej 

En 



T = 



En 



(48) 



(49) 



Let us consider the relation between the dimensionless 
temperature T and the absolute temperature T for the 
simple example of vitreous silica, the amorphous state of 
SiC>2 ■ Taking a microscopic length scale uq — lCP 10 m 
and substituting the values of h, ks and the mass m of 
SiC>2 in the definitions above, implies for this case that 
T = 1 corresponds approximately to T = 1 K. In this 
context we shall from here on look at the dimensionless 
temperature T as an approximate representation of the 
absolute temperature. 

In what follows we shall ignore writing the tildes on the 
dimensionless variables and parameters. One can show 
that the scaled form of all equations given in the previous 
sections is then obtained by setting h = m = 1. 



B. Matsubara correlations 

1. Perturbative solutions at two-loop order 

The perturbative Matsubara correlations were com- 
puted at two-loop order in the RS approximation. This 
requires solving the saddle-point equations f|32[) with 
qaa(T — t') = q c i(T — T ') which is an ft/3-periodic function, 
and q a ^b(T — T') = q. We treated them in a Fourier trans- 
formed representation. Our convention for the Fourier 
transform of fi./3-periodic functions /(r) is as follows 



(50a) 



f(u k ) = ^-JdT e-^W/(r), (50b) 

with Matsubara frequencies lq-^ — j^k (k — 0, ±1, ..). 
The advantage of this convention is that the dimension 
of the transformed quantity is equal to its original dimen- 
sion. The Fourier transform of the saddle-point equations 
(|32|) requires computing the Fourier transform q~ b (u>k) of 
the (functional) inverse kernel q~^(T, t'). Note that the 
latter defines an inverse w.r.t. both replica structure and 
Matsubara-time integration satisfying 



E 



/ dr"q- c }(r 



r")q cb (T"-r')=6 ab S(T-r'). (51) 
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Using Fourier transform relations requires that 

C'Wtot^) = 5 ab /(hf3) 2 , (52) 

c 

for all k, i.e. up to a factor \/(h(3) 2 the Fourier trans- 
forms <7ac(wfc) of the inverse kernel are equal to the cor- 
responding elements of the matrix inverse q(u>k)ac of the 
q{u>k) matrix in replica space, 



q- b \u k )=q(u Jk )- b 1 /m 2 



(53) 



The RS representation of the matrix elements q(ujk) a c m 
the n — > limit are given by [l8| 



1 



g(wfc) 



9d(wfe) - g(wfc) (qd(uk) - q(^k)) 2 
q{uk) 



(qd(u k ) -?(w fc )) 2 ' 



, (54a) 
(54b) 



where q(u>k) — qSk,o- Using (|5"4")) and |53^) in the Fourier 
transform of ((32)) leads to the following equations for the 
two-loop RS Matsubara correlations. First, for k — we 
have of 



1 



qd(uo) - q (qd(uo) - q) 2 



= - 2(/3J) 2 q d ^o) 

+ 24 5 /3^g d (w fe ), (55a) 



(?d( w o) - qf 



2(pJ) 2 q, 



(55b) 



where (|55a[) represents the replica-diagonal and (|55b[) the 
replica off-diagonal case. For /: / we only have to 
consider the replica-diagonal case, giving 



: pu 2 - 2(l3J) 2 q d {u Jk ) + 24g0j2^k) (56) 



qdi^k) 



In the high-temperature phase where q = the qd{<^k) 
are found numerically from (|55[) and (156|) by solving 
a single self-consistency equation for the variable i\ — 
J2k<ld(oJk) = <7d(0), namely 



(57) 



in which Bi. 



'2Agz\. In the low-temperature phase 



where q ^ 0, i\ can be determined analytically, entailing 
that the qd(^k) can be expressed in closed form as 



i 



(58a) 



(a**0 =4^2 (Bfc - ^B 2 -%.P) , (58b) 

(58c) 



V2 J 
1 12 s ' 




FIG. 2: Glass transition temperatures T s vs. J 1 at two-loop 
perturbation theory for g — 1. 



The glass transition temperature T g as a function of 
J (at g = 1) is shown in Fig. [51 Glass transition tem- 
peratures for structural glasses are typically in the range 
between 500 K and 1500 K. Clearly this requires J to be 
large. The majority of our results in the present study 
were therefore computed for a typical large J, J = 50 
giving T g w 500 K. 

Next is to discuss the solutions of (|58[) for the Mat- 
subara correlations. The gd(r) were computed from a 
numerical inverse Fourier transformation of qdfak)- In 
Fig. [3] we plot the results for qd(i~) — g for a number of 
low temperatures. A selection of them is compared with 
the results of Quantum Monte Carlo (QMC) simulations 
of the ncm-perturbative theory, which will be discussed 
in a subsection below. The solutions for the Matsubara 




FIG. 3: Two-loop perturbative data for the function qd(r) — 
q. Inverse temperatures (in K^ 1 ) from top to bottom are 
/3 = 0.1, 0.13, 0.2, 0.3, 0.4, 0.5, 1, 2, 5, 10, 50. Comparison with 
QMC solutions of the non-perturbative theory for (3 = 0.5 
and /3 = 1 (marked by + and x). 

correlation q are plotted in Fig. |4]in the low-temperature 
phase. From both Fig. [3]and[5]we conclude that the two- 
loop results are in reasonably good agreement with the 
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FIG. 4: Perturbative and non-perturbative data for the Mat- 
subara correlation q in the low-temperature phase. Note 
the good agreement of two-loop perturbative data with non- 
perturbative data from Quantum Monte Carlo simulations. 
Data of non-perturbative static approximations deviate from 
these in the low-temperature region. 



QMC results. 

The Fourier-representation of the self-consistency 
equations for the Matsubara correlations in the low- 
temperature phase at three-loop order is given in the 
appendix. The coupling between the Fourier modes at 
three-loop order contains a truly functional element via 
the function 23 (u>k)- Although we succeeded in simpli- 
fying the problem to solving a set of only four coupled 
transcendental equations for the variables 2i,^2>9d( w o) 
and q, we have so far been unable to solve them. In fact 
we suspect that physically acceptable solutions may not 
exist at three-loop order and an expansion to higher loop 
order might be necessary. 



2. Non-perturbative Quantum Monte Carlo simulations 

The Matsubara correlations <7d(r) and q of the non- 
perturbative RS theory defined in section IVI Al above, 
were evaluated with Quantum Monte Carlo simulations. 
This involved solving the functional self-consistent rela- 
tions (I40|). We used iterative QMC-techniques along the 
lines of [23j . starting with a set of initial values of qd(T~) 
and q to be used as input for the action ([57]) . after which 
they were updated in a path integral Monte Carlo algo- 
rithm. This procedure was repeated 10 times, resulting 
in reasonably good convergence of the Matsubara correla- 
tions. As regards to the algorithm, an update contained 
10 5 Monte Carlo sweeps (taking data every 10th sweep), 
5 • 10 4 equilibration sweeps and 5 • 10 3 Gaussian z, z sam- 
ples. The imaginary time axis was discretised into 40 
time slices. The results for — q at a selected num- 

ber of temperatures are plotted in Fig. [3l The results 
for the off-diagonal Matsubara correlation q are plotted 
separately in Fig. [5] As mentioned previously, they were 
found to be in good agreement with the solutions of the 



two- loop perturbative theory given in (|58[) . The conclu- 
sion from the simulations is that they confirm the validity 
of the perturbative two-loop results for the Matsubara 
correlations. 



3. Non-perturbative static solutions 

The Matsubara correlations of the static approxima- 
tion treated in section fVI B[ were computed numerically. 
This involved solving the functional self-consistency re- 
lations (|45[) and (|46[) . They were solved in the opera- 
tor representation, for which the required Hamiltonian 
H s t(p,u,v;z,z) is reconstructed from the action (|4"2j) . 
Path integrals are then re-expressed in terms of traces 
over a suitable truncated Hilbert space. We used a 
bases of harmonic oscillator eigenstates. The Gaussian 
integrals integrals were computed with Gauss-Legendre 
quadratures. 

The results for the replica off-diagonal Matsubara cor- 
relations q are plotted in Fig. Uover a large temperature 
range. We believe the differences with the two-loop per- 
turbative and QMC results seen here to be an artifact of 
the static approximation, which after all does not repre- 
sent the true self-consistent solutions of the saddle-point 
equations (|4TJ)l . 

Differences between the results of the two- variable and 
three-variable approach are not visible on the scale used 
in Fig. [4] They are plotted separately in Fig. [5] for a 
selected range of low temperatures. The combination of 
order parameters C — (3(qd — q), which can be seen as 
a susceptibility-like variable, was found to be approxi- 
mately equal for both, the two-variable and the three- 
variable approach, with nearly constant numerical value 
C re 17 • 10~ 3 for < T < 100 K. Since this is very 
small, the values of qd and q barely differ at low temper- 
atures, as is indicated in Fig. O On the other hand, the 
two-variable and three-variable approaches do give rise 
to different values for the (qd, q) pairs: Introduction of a 
third variable ^(O) in the three- variable approach leads 
to a depression of the values of qd and q relative to those 
in the two- variable approach, whereas <?d(0) in the three- 
variable approach turns out to be larger than qd and q 
within the two-variable solution. 



C. Thermodynamics 

1. Perturbative specific heat 

To obtain an expression for the perturbative free en- 
ergy we substitute (frj| and in (|5U]). We consider 
Fourier transforms as defined in (|50p and the scaling in- 
troduced in section [VII Al The trace ^Trg^g diverges 
when substituting q^ 1 from its definition (| 1 T|) . Instead 
we substitute the expression for q^ 1 determined by the 
saddle-point equations ([52)1 or ([54"]) . The result is then 
finite up to an irrelevant infinite constant ^Trq^ 1 ^. The 
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3.96 



0.01 



T(K) 



FIG. 5: The upper line shows qd(0) in the three- variable ap- 
proach. The lowest pair of lines correspond to qd and q in 
the three- variable approach, whereas the pair of lines in the 
middle correspond to qd and q in the two- variable approach. 



trace iTrlogq 1 can be shown to have the following RS 



representation [18| 



o 
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FIG. 6: Approx. linear specific heat C ~ T 102 (T < 0.8 iT) 
and C ~ T 1 1 (T < 0.5 if) for resp. perturbative (three- 
loop) and non-perturbative (static) theory and a Bose-peak 
at higher temperatures. The two-loop perturbative theory 
described only the vibrational excitations. 



2. Non-perturbative specific heat 



-Ty loe j -1 = "Vf ^ fc ) 

2 1 gq 2^\q d {^k)-q(0Jk) 

+ \og(q d (uj k )-q(LJ k ))}, (59) 



where n defines the number of replicas. Observe that 
when substituting q{uj k ) — qd k .o the second part of (j59"|) 
contains the divergent sum — | Ylkj^o l°S9rf( w fe)' To deal 
with this we first substitute gd(wfe) from (|56[) or (p2^) . 
From the result we isolate a divergent contribution of 
the following form ^J2 k log {(3(^1 + 24gzi)} with zj = 
Sfe Qd(wk)- This term, as part of (|59|) and r o ff in ([28]) . 
should be exponentiated according to (p?Tj) , defining the 
partition function of a simple harmonic oscillator with 
frequency ojq = y/2Agz\. Consequently, we may replace 
\ J2k 1°S {P&k + w o)} by (/3 times) the free energy of a 
simple harmonic oscillator, giving logsinh(i/3w ) which 
is now finite. 

To compute the RS free energy numerically we have 
only the two-loop data from ([55)1 at our disposal. The 
specific heat computed from the two-loop free energy did 
not result in a glassy low-temperature anomaly. Only 
vibrational excitations of the system featured here (see 
Fig. [5]). On the other hand, when evaluating the free 
energy at three-loop order, a specific heat exhibiting the 
characteristic glassy low-temperature anomaly was ob- 
tained (C ~ T lm for T < 0.8 K), though we had to use 
two-loop results for the Matsubara correlations in those 
expressions (as three-loop results are so far unavailable). 
The good agreement between two-loop Matsubara corre- 
lations and QMC results is thought to provide a reason- 
able justification for this approach. 



The non-perturbative thermodynamics was evaluated 
in the static approximation. For this we used the free 
energy expression (|44| and the numerical results de- 
termined from (|45]) and (|46|) . This indeed reproduced 
the characteristic glassy low-temperature specific heat 
anomaly for both the two- variable and three- variable ap- 
proach (see Fig. [6|). Again the low-temperature spe- 
cific heat showed an approximately linear (in fact super- 
linear) temperature dependence C ~ T 11 for T < 0.5^, 
in reasonable agreement also with experimental data 
[Uliij]. These results should also be compared with, and 
are indeed comparable to those of the translationally in- 
variant model investigated in Q. We found little differ- 
ence in the results for the two- variable and three- variable 
approach at higher temperatures, as can be seen in Fig. 
[6] Differences between the three-loop perturbative spe- 
cific heat and the non-perturbative 'static' specific heat 
are restricted to the 0.5 — 5 K temperature region. This 
could be at least partly because of the strong tempera- 
ture dependence of the 'static' Matsubara correlations in 
this region, as displayed in Fig. |4]for the order parameter 

q- 

Finally we remark that both the properties of the Bose- 
peak at intermediate temperatures and the universal tun- 
neling regime at low temperatures are governed by one 
and the same set of system parameters appearing in (|42[) . 
No separate sets of assumptions were introduced to de- 
scribe these two temperature regimes. 



VIII. CONCLUSIONS 

In summary, we have provided a fully quantum statisti- 
cal analysis of a microscopic model of a glass, respecting 
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global translation invariance. Until now such analysis 
was available only at a semi-classical level. We formu- 
lated an effective theory in terms of single-site path inte- 
grals and constructed perturbative and non-perturbativc 
solutions of a set of self-consistency equations describing 
the system. Both resulted in an approximately linear spe- 
cific heat at low temperatures, in good agreement with 
experiment. 

The perturbation theory was formulated in terms of 
two-particle irreducible diagrams at two-loop and three- 
loop order for the effective action. As solutions of the 
self-consistency equations at three-loop order remained 
unavaible, we resorted to investigating the reliability of 
our two-loop results using Quantum Monte Carlo sim- 
ulations. We found surprisingly good agreement of the 
Matsubara correlations obtained perturbatively and via 
QMC simulations. 

Within a non-perturbative static approximation we ob- 
tain a description in terms of a glassy potential energy 
landscape containing an ensemble of effective single-well 
and double-well potentials, much as in the soft-potential 
model 0, H| and in the semi-classical approach In- 
terestingly there is an important difference, namely the 
emergence of a coupling to an additional classical variable 
in a manner reminiscent of a coupling of local excitations 
to a heat bath as postulated within phenomenological 
models 

It would be interesting to carry the perturbative ap- 
proach to higher loop order, or in fact attempt summa- 
tions of infinite classes of 2PI diagrams. On another 
front, effects of replica symmetry breaking have not yet 
been looked at and are worth investigating (though in 
a semi-classical approach RSB effects were found to be 
weak [§])• 

One of the motivations for the present investigation 
was to understand a phase transition observed in ultra- 
cold glasses more than a decade ago (2(| , which has so- far 
not found an explanation. Regrettably, the present study 
has not produced any progress in that particular direc- 
tion. It might well be the case that an expansion of the 
present investigation in both directions mentioned above 
— including effects of replica symmetry breaking and in- 
clusion of diagrams up to arbitrarily high loop order — 
would be required to reveal pertinent signatures of that 
phase transition. 



IX. APPENDIX 

Here we present the Fourier transformed RS represen- 
tation of the three-loop perturbative equations (134|) in 
the low-temperature phase. Considering the scaling in- 
troduced in section IVII Al and using ([55]) and ([51)1 , the 
Fourier transformation leads to the following set of equa- 
tions. First, for k = we have 

-r-— r - - ^ = - 2{[3J) 2 q d {u ) + 2Agpz 1 

qd{uo)-q (qd(^o)-qr 

-96((3g) 2 z 3 (uj ) 

+ 24(f3J) 2 gq d (LO )f3(z 2 -q 2 ), 
(60a) 

-T-tA — ^ = - 2(/3J) 2 g 
(q d {uJo) - q) 2 

-96(Pg) 2 q 3 

+ 24(pj) 2 gqP(z 2 -q 2 ), 
(60b) 

where (|60a[) represents the replica-diagonal and (|60b[) the 
replica off-diagonal case. The last two terms in these 
equations define the three-loop extensions of the two- loop 
equations (|55|) . Here we used the following definitions 

z 2 =^qd(^k) 2 , (61a) 

k 

h{uk) = <ld(ui) qd(u m ) (jdjuk -t*>j -w m ). (61b) 

im 

For k ^ we only have to consider the replica-diagonal 
case, giving 

— l — = Pu 2 -2(pj) 2 q d (uk) + 24gP Zl 
qd{cok) 

-96((3g) 2 z 3 (cj k ) 

+24(pj) 2 gq d (uk)P(z 2 ~q 2 ), (62) 

where the last two terms again represent the three-loop 
extensions of the two-loop equations. The solutions of 
(I60p can no more be expressed in analytic form as was 
the case for the two-loop perturbative equations ([SS]) . 
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